#
# Replication file for 
# "Leave It As It Is": International Network Effects on Protected Lands
# Darren Hawkins and Jay Goodliffe
# International Interactions
#
# Social Network Analysis using Siena 
#

library(RSiena)

# steps: 
# export protected lands and domestic variables as separate wide Stata files 
# import into R, make into martrix file, and into Siena files
# network variable: create edge matrices using cutoffs (country x country for each year)
#                   (OK as Stata files)
# import into R and use igraph to make graph file, then adjacency files 
#                   (for each year, each network)
# make into sienaDependent variables stacking the years

setwd("directory path here")
load("II Replication Data.RData")

# Trade Network
pl_data <- sienaDataCreate(pl3.siena, trade1.siena, polity.siena, gdppc.siena, 
                           natrent.siena, popden.siena, ruralpop.siena, rulelaw.siena)
pl_data
pl_eff <- getEffects( pl_data )
# standard network effects
pl_eff <- includeEffects( pl_eff, name = "trade1.siena", transTrip, inPop )
# effects of network on protected lands
pl_eff <- includeEffects(pl_eff, name = "pl3.siena", avSim, interaction1="trade1.siena")
# drop shape parameters, as no global (collinear) effect
pl_eff <- setEffect(pl_eff, linear, include = FALSE, name = "pl3.siena", type = "eval")
pl_eff <- setEffect(pl_eff, quad, include = FALSE, name = "pl3.siena", type = "eval")
# covariate effects
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "polity.siena" )
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "gdppc.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "natrent.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "popden.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "ruralpop.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "rulelaw.siena")
# effect of similarity in protected lands on network formation
pl_eff <- includeEffects(pl_eff, name = "trade1.siena", simX, interaction1="pl3.siena")
pl_eff
# set up algorithm
pl_algorithm <- sienaAlgorithmCreate( projname = 'protectedlands' )
# Estimate parameters
pl_results_trade <- siena07( pl_algorithm, 
                                  data = pl_data, 
                                  effects = pl_eff, 
                                  useCluster = TRUE, #if using parallel processing
                                  nbrNodes = 9, #if using parallel processing
                                  initC = TRUE)
pl_results_trade

# if have not yet met convergence, continue
siena07ToConvergence <- function(alg, dat, eff) {
  numr <- 1
  ans_t <- pl_results_pl3_trade1 # the previous run
  repeat {
    save(ans_t, file=paste("ans_t_",numr,".RData",sep="")) # to be safe
    numr <- numr+1 # count number of repeated runs
    tm <- ans_t$tconv.max # convergence indicator
    cat(numr, tm,"\n") # report how far we are
    if (tm < 0.25) {break} # success
    if (tm > 10) {break} # divergence without much hope
    # of returning to good parameter values
    if (numr > 30) {break} # now it has lasted too long
    ans_t <- siena07(alg, 
                     data=dat, 
                     effects=eff, 
                     useCluster = TRUE, 
                     nbrNodes = 9, 
                     initC = TRUE,
                     prevAns=ans_t)
    ans_t
    numr
  }
  if (tm > 0.25)
  {
    cat("Warning: convergence inadequate.\n")
  }
  ans_t
}
siena07ToConvergence(pl_algorithm, pl_data, pl_eff)

# replace N with last convergence (7?)
load("ans_t_N.RData")
pl_results_trade_final <- ans_t


# IO Network
pl_data <- sienaDataCreate(pl3.siena, igo10_6.siena, polity.siena, gdppc.siena, 
                           natrent.siena, popden.siena, ruralpop.siena, rulelaw.siena)
pl_data
pl_eff <- getEffects( pl_data )
# standard network effects
pl_eff <- includeEffects( pl_eff, name = "igo10_6.siena", transTrip, inPop )
# effects of network on protected lands
pl_eff <- includeEffects(pl_eff, name = "pl3.siena", avSim, interaction1="igo10_6.siena")
# drop shape parameters, as no global (collinear) effect
pl_eff <- setEffect(pl_eff, linear, include = FALSE, name = "pl3.siena", type = "eval")
pl_eff <- setEffect(pl_eff, quad, include = FALSE, name = "pl3.siena", type = "eval")
# covariate effects
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "polity.siena" )
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "gdppc.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "natrent.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "popden.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "ruralpop.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "rulelaw.siena")
# effect of similarity in protected lands on network formation
pl_eff <- includeEffects(pl_eff, name = "igo10_6.siena", simX, interaction1="pl3.siena")
pl_eff
# set up algorithm
pl_algorithm <- sienaAlgorithmCreate( projname = 'protectedlands' )

# Estimate parameters
pl_results_io <- siena07( pl_algorithm,
                           data = pl_data, 
                           effects = pl_eff, 
                           useCluster = TRUE, #if using parallel processing
                           nbrNodes = 9, #if using parallel processing
                           initC = TRUE)
pl_results_io

# if have not yet met convergence, continue
siena07ToConvergence <- function(alg, dat, eff) {
  numr <- 1
  ans_i <- pl_results_io # the previous run
  repeat {
    save(ans_i, file=paste("ans_i_",numr,".RData",sep="")) # to be safe
    tm <- ans_i$tconv.max # convergence indicator
    cat(numr, tm,"\n") # report how far we are
    numr <- numr+1 # count number of repeated runs
    if (tm < 0.25) {break} # success
    if (tm > 10) {break} # divergence without much hope
    # of returning to good parameter values
    if (numr > 30) {break} # now it has lasted too long
    ans_i <- siena07(alg, data=dat, effects=eff, 
                     useCluster = TRUE, nbrNodes = 9, initC = TRUE,
                     prevAns=ans_i)
  }
  if (tm > 0.25)
  {
    cat("Warning: convergence inadequate.\n")
  }
  ans_i
}
siena07ToConvergence(pl_algorithm, pl_data, pl_eff)

# replace N with last convergence (3?)
load("ans_i_N.RData")
pl_results_io_final <- ans_i


# Aid Network
pl_data <- sienaDataCreate(pl3.siena, aid05.siena, polity.siena, gdppc.siena, 
                           natrent.siena, popden.siena, ruralpop.siena, rulelaw.siena)
pl_data
pl_eff <- getEffects( pl_data )
# standard network effects, but no transTrip
pl_eff <- includeEffects( pl_eff, name = "aid05.siena", inPop )
# remove aid reciprocity
pl_eff <- setEffect(pl_eff, recip, include = FALSE, name = "aid05.siena", type = "eval")
# set the missing rate for aid
pl_eff <- setEffect(pl_eff, Rate, fix = TRUE, initialValue = 10, name = "aid05.siena",
                    type = "rate", period = 26)
# effects of network on protected lands
pl_eff <- includeEffects(pl_eff, name = "pl3.siena",avSim, interaction1="aid05.siena")
# drop shape parameters, as no global (collinear) effect
pl_eff <- setEffect(pl_eff, linear, include = FALSE, name = "pl3.siena", type = "eval")
pl_eff <- setEffect(pl_eff, quad, include = FALSE, name = "pl3.siena", type = "eval")
# covariate effects
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "polity.siena" )
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "gdppc.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "natrent.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "popden.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "ruralpop.siena")
pl_eff <- includeEffects(pl_eff, effFrom, name = "pl3.siena", interaction1 = "rulelaw.siena")
# effect of similarity in protected lands on network formation
pl_eff <- includeEffects(pl_eff, name = "aid05.siena", simX, interaction1="pl3.siena")
pl_eff
# set up algorithm
pl_algorithm <- sienaAlgorithmCreate( projname = 'protectedlands' )

# Estimate parameters
pl_results_aid <- siena07( pl_algorithm,
                          data = pl_data, 
                          effects = pl_eff, 
                          useCluster = TRUE, #if using parallel processing
                          nbrNodes = 9, #if using parallel processing
                          initC = TRUE)
pl_results_aid

siena07ToConvergence <- function(alg, dat, eff) {
  numr <- 1
  ans_a <- pl_results_aid # the previous run
  repeat {
    save(ans_a, file=paste("ans_a_",numr,".RData",sep="")) # to be safe
    tm <- ans_a$tconv.max # convergence indicator
    cat(numr, tm,"\n") # report how far we are
    numr <- numr+1 # count number of repeated runs
    if (tm < 0.25) {break} # success
    if (tm > 10) {break} # divergence without much hope
    # of returning to good parameter values
    if (numr > 30) {break} # now it has lasted too long
    ans_a <- siena07(alg, data=dat, effects=eff, 
                     useCluster = TRUE, nbrNodes = 9, initC = TRUE,
                     prevAns=ans_a)
  }
  if (tm > 0.25)
  {
    cat("Warning: convergence inadequate.\n")
  }
  ans_a
}
siena07ToConvergence(pl_algorithm, pl_data, pl_eff)

# replace N with last convergence (24?)
load("ans_a_24.RData")
pl_results_aid_final <- ans_a

save.image(file = "protectedlands_eachnet.RData")
